#include "testmv.h"
#include "em_fit.h"

/*
THESE TESTS DO NOT CHECK ANYTHING!

They are NOT tests, just experimentation. This whole file should be reworked.
*/

/*
n=300;
x = sqrt(randn(n)**2 + randn(n)**2);
hist(x,bins=2*sqrt(n));
1./2./n*sum(x**2) # estimated variance

m=300;
y = rand(m,2); # rejection sampling for triangle
y = y[y[:,1] < y[:,0], 0] * 5

hold(True)
hist(y,bins=2*sqrt(m),fc='r');

hstack((x,y)) # generate below data
*/

double test_values[] = {
        0.56521647,  0.33243875,  1.28457487,  1.32644918,  1.40217558,
        1.72814282,  0.29182806,  1.58979241,  1.31693304,  1.27134406,
        0.69182237,  2.30068956,  1.12119675,  0.8040021 ,  2.92711315,
        1.28704531,  0.45112437,  0.59897044,  0.92212769,  0.65613286,
        1.69755479,  1.64886821,  0.57443999,  1.59149849,  0.63194567,
        0.23519801,  1.42118414,  0.47934526,  1.5618863 ,  1.16268465,
        1.42559416,  1.06520252,  0.83384824,  0.94545045,  1.63084661,
        4.37052128,  0.69893292,  0.61506276,  0.81580177,  2.24706533,
        1.45086483,  1.12135264,  0.84729055,  1.41254275,  0.89751047,
        1.81606689,  2.24905105,  1.12125889,  0.67771351,  2.37280699,
        1.91404933,  0.74139973,  2.23036143,  2.38619772,  1.81922661,
        1.38698507,  0.88929874,  0.79695151,  1.39879063,  0.48311405,
        0.70184236,  1.52672943,  2.17922692,  0.79925609,  1.69775721,
        1.79605353,  0.95416743,  0.49022438,  0.79918622,  1.61786289,
        2.80158992,  2.02143135,  1.86006584,  0.96511829,  1.74495942,
        0.7117868 ,  1.82816142,  0.66774918,  0.55613942,  0.72658459,
        0.50142717,  0.76915237,  1.72635425,  1.08177112,  0.24652656,
        0.5211352 ,  1.24130487,  1.93571071,  0.4365494 ,  2.10912617,
        0.30681448,  2.87463332,  1.33949035,  0.654821  ,  1.71437981,
        1.4693267 ,  0.78775527,  0.08177178,  0.32417893,  1.6994855 ,
        0.38586525,  1.21946215,  0.87850911,  1.28873245,  0.72654688,
        0.70628648,  1.29998193,  1.43081093,  0.68249697,  1.71770923,
        0.72752214,  1.80593291,  1.31492758,  0.61569317,  1.92535144,
        1.00946003,  1.97702079,  0.3105334 ,  0.69733709,  0.59327191,
        1.47846489,  1.01696555,  0.9983398 ,  1.46881122,  2.0947015 ,
        1.23525646,  2.14541564,  0.59770355,  0.52778717,  1.00606361,
        1.65786651,  0.75531171,  1.4097308 ,  0.72253515,  1.81774273,
        0.43880844,  0.62625475,  0.61746296,  0.16988884,  0.89170407,
        0.81077139,  1.51833464,  0.87132356,  1.33793503,  1.41008224,
        0.86634915,  0.94758522,  0.65224458,  1.95046595,  1.28555969,
        1.42214334,  0.72859512,  0.66181897,  1.16975019,  0.62254101,
        1.69163761,  0.89149653,  1.42762242,  1.53249792,  2.39245171,
        0.55645568,  1.02760322,  1.39730257,  0.69637063,  0.51315557,
        1.18445672,  1.12205516,  1.67322178,  0.80039271,  1.39203534,
        0.44336421,  3.821907  ,  1.25926231,  0.65036715,  2.97116038,
        0.82970407,  0.1825084 ,  1.28972257,  0.50582226,  1.34462737,
        1.41088449,  0.95509973,  1.06921248,  0.51741079,  2.45501073,
        1.60281862,  1.04867476,  0.8398178 ,  1.8984844 ,  1.10188191,
        2.05793156,  0.09287777,  0.687779  ,  1.33442868,  0.50280667,
        2.10685789,  1.96116359,  2.2685843 ,  2.462193  ,  0.72898301,
        1.02616034,  1.44343835,  0.80208236,  0.69412428,  1.13324193,
        0.92419159,  1.46275601,  0.89706023,  2.36033676,  0.70574998,
        1.28256802,  1.08997899,  1.1537689 ,  1.50893188,  0.97297716,
        1.82313152,  0.08819429,  1.91391277,  0.35388534,  1.43213824,
        2.4382917 ,  1.95027624,  0.76685544,  2.40667373,  1.2812412 ,
        2.43411992,  1.47760654,  1.6749107 ,  1.53207962,  1.34379747,
        0.89138885,  0.21513171,  2.21322507,  1.9865672 ,  1.02324029,
        2.2415141 ,  0.35573617,  1.77589608,  0.63552525,  2.16300654,
        0.86327281,  0.62409708,  1.75655824,  2.53665564,  1.73830601,
        0.60320183,  2.30381739,  0.55724155,  0.69195099,  1.52449512,
        0.52992514,  1.29542693,  1.5928496 ,  0.9240627 ,  1.23248694,
        1.15222835,  1.25984132,  1.29683717,  1.66067196,  1.90634925,
        1.51810124,  1.11765724,  0.84156329,  1.95480656,  1.10715653,
        1.17168321,  0.98470028,  0.6050028 ,  1.74792011,  1.13158615,
        0.46581404,  0.87708388,  0.69477858,  0.72407576,  2.08172137,
        1.10177336,  1.73029337,  1.62842499,  1.2597808 ,  1.94174256,
        1.25309353,  0.71696458,  1.63417676,  0.76105528,  0.78794606,
        2.45597183,  1.37655407,  0.44743302,  0.95965156,  2.20776749,
        0.70998322,  1.5027764 ,  1.68259952,  0.63991464,  0.54674961,
        0.639517  ,  0.78233136,  1.20898756,  1.5895498 ,  1.38142862,
        3.59977092,  2.71742438,  0.2638367 ,  3.21456485,  2.27763291,
        3.2702455 ,  4.82999047,  3.97495639,  4.92647734,  4.30317769,
        4.19251638,  4.79527369,  3.35389452,  4.12874021,  2.72369012,
        3.55756225,  2.18618109,  2.08466423,  3.62988934,  3.51855219,
        0.72115977,  4.71015019,  2.53459855,  2.92627837,  4.63452221,
        1.57302468,  4.31611432,  3.86302858,  3.0753069 ,  3.0857125 ,
        3.13916188,  4.46148803,  3.40529874,  3.18041005,  3.76428837,
        3.3640308 ,  0.35603387,  3.9513179 ,  2.54107137,  1.48793218,
        3.81370287,  4.3919535 ,  4.80873303,  3.47925805,  2.97557602,
        2.53947533,  2.60510082,  1.96939411,  4.95334435,  3.50664232,
        2.15066291,  3.66983785,  4.53939282,  2.35490479,  4.49764476,
        2.97865863,  2.45814149,  4.87742553,  2.89082454,  3.79711842,
        3.47145118,  3.86658496,  3.26087407,  4.83164916,  3.62180073,
        0.84772602,  4.60487853,  4.96135533,  3.93109223,  4.23854872,
        3.11217977,  4.2211854 ,  4.69793729,  0.83517887,  4.29903205,
        4.51170977,  3.95408048,  4.41162494,  4.64169263,  1.75649843,
        3.11480628,  4.0172625 ,  2.60356951,  2.64617035,  2.52016811,
        2.12705187,  1.46959111,  2.6345852 ,  3.92159828,  2.80085297,
        2.30835482,  3.14298989,  3.16268608,  2.9655581 ,  2.31281831,
        4.99038274,  2.61541704,  4.17467769,  3.41851594,  4.87184856,
        4.28298523,  4.49014588,  0.84513533,  4.14465731,  4.04439945,
        4.42856269,  4.8499938 ,  1.74225175,  0.87595661,  4.13856459,
        3.45701148,  2.86594465,  2.12907041,  3.49183867,  4.40023839,
        4.66458584,  2.5328817 ,  2.82132741,  4.23112263,  2.22985953,
        2.77874983,  4.31329982,  3.33386324,  4.12260939,  4.69871697,
        1.5212395 ,  2.44277537,  3.98636126,  4.80118806,  3.24456176,
        3.06643682,  0.88130849,  2.5763712 ,  2.82121838,  3.88340189,
        2.83671598,  2.20218383,  3.43659774,  4.60316828,  4.7797459 ,
        2.42694069,  2.92767903,  4.01287253,  2.29269741,  4.77270134,
        2.13799981,  4.93933666,  3.7014534 ,  4.80479814,  2.17039502,
        1.95943342,  3.84809163,  4.34025879,  4.01615128,  1.46365492,
        4.16528947,  2.53096174,  3.254525  ,  1.64440713,  3.93587379
};

TEST() {
	vec x(2);
	x[0] = 10.5;
	x[1] = 3.0;
	Nearby(logsum(x), 10.5005529315, 1e-6);
}

TEST() {
	RayleighDistribution r(0.2);
	std::vector<Distribution*> distribs;
	distribs.push_back(&r);

	int nsamps = sizeof(test_values)/sizeof(double);
	// first 300 samples are pure rayleigh
	nsamps = 300;
	vec residuals(nsamps);
	for (int i=0; i<nsamps; i++)
		residuals[i] = test_values[i];

	vec drops = vzeros(nsamps);
	double lp = em_fit(residuals, distribs, drops);
	printf("lp=%g\n",lp);
	printf("sigma2 r1 %g\n",r.sigma2);

	Nearby(r.sigma2, 1.0, 0.1);
}

TEST() {
	RayleighDistribution r(0.1);
	RightUniformDistribution u(20);
//	RayleighDistribution u(10.2);
	std::vector<Distribution*> distribs;
	distribs.push_back(&r);
	distribs.push_back(&u);

	int nsamps = sizeof(test_values)/sizeof(double);
//	nsamps = 300;
	vec residuals(nsamps);
	for (int i=0; i<nsamps; i++)
		residuals[i] = test_values[i];

	vec drops = vzeros(nsamps);
	double lp = em_fit(residuals, distribs, drops);
	printf("sigma2 r1 %g\n",r.sigma2);
//	printf("sigma2 r2 %g\n",u.sigma2);
	printf("a %g\n",u.a);
	printf("lp=%g\n",lp);
}

TEST() {
	RayleighDistribution r(0.1);
//	RightUniformDistribution u(20);
	RayleighDistribution u(10.2);
	std::vector<Distribution*> distribs;
	distribs.push_back(&r);
	distribs.push_back(&u);

	int nsamps = sizeof(test_values)/sizeof(double);
//	nsamps = 300;
	vec residuals(nsamps);
	for (int i=0; i<nsamps; i++)
		residuals[i] = test_values[i];

	vec drops = vzeros(nsamps);
	double lp = em_fit(residuals, distribs, drops);
	printf("sigma2 r1 %g\n",r.sigma2);
	printf("sigma2 r2 %g\n",u.sigma2);
//	printf("a %g\n",u.a);
	printf("lp=%g\n",lp);
}

TEST() {
	printf("////////////////////////////////////////////\n");
	RayleighDistribution r1(0.1);
	RayleighDistribution r2(10.2);
	RightUniformDistribution u(20);
	std::vector<Distribution*> distribs;
	distribs.push_back(&r1);
	distribs.push_back(&r2);
	distribs.push_back(&u);

	int nsamps = sizeof(test_values)/sizeof(double);
//	nsamps = 300;
	vec residuals(nsamps);
	for (int i=0; i<nsamps; i++)
		residuals[i] = test_values[i];

	vec drops = vzeros(nsamps);
	double lp = em_fit(residuals, distribs, drops);
	printf("sigma2 r1 %g\n",r1.sigma2);
	printf("sigma2 r2 %g\n",r2.sigma2);
	printf("a %g\n",u.a);
	printf("lp=%g\n",lp);
}


TEST() {
	printf("////////////////////////////////////////////\n");
	RightUniformDistribution u1(100);
	RightUniformDistribution u2(20);
	std::vector<Distribution*> distribs;
	distribs.push_back(&u1);
	distribs.push_back(&u2);

	int nsamps = sizeof(test_values)/sizeof(double);
//	nsamps = 300;
	vec residuals(nsamps);
	for (int i=0; i<nsamps; i++)
		residuals[i] = test_values[i];

	vec drops = vzeros(nsamps);
	double lp = em_fit(residuals, distribs, drops);
	printf("a %g\n",u1.a);
	printf("a %g\n",u2.a);
	printf("lp=%g\n",lp);
}

double tricky_residuals[] = {3.78891e-08, 2.96607e-08,
	3.57762e-08, 4.20602e-08, 9.55581e-10, 1.21278e-07,
	9.59277e-08, 5.75477e-08, 8.9757e-08, 2.483e-08, 2.62061e-08,
	9.13176e-09, 9.46357e-09, 6.19956e-10, 2.58106e-08,
	1.60917e-08, 4.17096e-08, 2.05079e-08, 5.74833e-08,
	2.73954e-09, 2.21134e-08, 5.28591e-08, 9.2879e-10,
	1.03962e-07, 1.38619e-07, 2.29725e-10, 1.00943e-07,
	1.9946e-08, 1.16861e-07, 2.71067e-10};

TEST() {
	printf("////////////////////////////////////////////\n");
	RayleighDistribution r1(0.1);
	RayleighDistribution r2(10.2);
	std::vector<Distribution*> distribs;
	distribs.push_back(&r1);
	distribs.push_back(&r2);

	int nsamps = sizeof(tricky_residuals)/sizeof(double);
	vec residuals(nsamps);
	fill(residuals, tricky_residuals);

	vec drops = vzeros(nsamps);
	double lp = em_fit(residuals, distribs, drops);
	printf("a %g\n",r1.sigma2);
	printf("a %g\n",r2.sigma2);
	printf("lp=%g\n",lp);
	printf("dropped=%g\n",sum(drops));
}

